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Abstract 

We discuss worm algorithms for the 3-state Potts model with ex- 
ternal field and chemical potential. The complex phase problem of 
this system can be overcome by using a flux representation where the 
new degrees of freedom are dimer and monomer variables. Working 
with this representation we discuss two different generalizations of the 
conventional Prokof'ev-Svistunov algorithm suitable for Monte Carlo 
simulations of the model at arbitrary chemical potential and evaluate 
their performance. 
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1 Introductory remarks 



Over the last two decades lattice methods have made considerable progress 
in our understanding of non-perturbative QCD, and now give rise to reliable 
quantitative results. However, one of the areas where progress is painfully 
slow is the analysis of QCD at finite density. The introduction of a chemical 
potential turns the fermion determinant into a complex number, such that 
it cannot be used as a weight factor in a Monte Carlo analysis - an obstacle 
known as complex phase problem. For true progress with QCD thermody- 
namics on the lattice new algorithmic concepts will be necessary. Often such 
new ideas are first developed and tested in simpler models before becoming 
applicable to the full theory. 

The worm algorithm [I] is a conceptually new approach which has al- 
ready been shown to be useful in some QCD related applications, in particu- 
lar in strongly coupled lattice QCD [21 [3] and in scalar [4] or low dimensional 
lattice field theories [5]. In this paper we extend the range of applications 
of worm algorithms to QCD related problems, and show that a worm algo- 
rithm can be used to simulate the 3-state Potts model in three dimensions 
with a center symmetry breaking term and a chemical potential, which in 
its standard representation has a complex phase problem. This model has 
been discussed as an effective theory for parts of the QCD phase diagram 
[HI [JJ [8] . Numerical studies with various techniques for vanishing [91 [10] and 
non- vanishing chemical potential [HI \TT\ [T2l [T3] may be found in the liter- 
ature. Recently a systematical analysis of the phase diagram as function 
of temperature and chemical potential was presented and its implications 
for QCD were discussed [14] . In this work we discuss in detail the worm 
algorithms used in [14] and assess their performance. 

The 3-state Potts model with magnetic field k and chemical potential [i 
in d dimensions is described by the Hamiltonian 

H[P] = ( T Y1 [P{x)P{x+i>)* + c.c] + rjP(x) + 7jP(x)*\ , (1) 

X \ u = l / 

where we use the variables r/ and 7] which are related to the strength k of 
the external field and the chemical potential [i via 

rj = Ke^ , rj = Ke~ M . (2) 

The first sum in (TTJ) is over all sites x of a d-dimensional (hyper-) cubic 
lattice with periodic boundary conditions and a total of V lattice points. 
The second sum is over all directions v = 1,2, ...d, where v denotes the 
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unit vector in direction v. The spin variables P (x) are complex phases 
with three different values, P(x) G Z3 = {l,e l2w ^ 3 ,e~ l27T ^ 3 }. The coupling 
parameters r, k and \i are real and non-negative. The partition function of 
the model is given by Z = J2{p} e~ S ^- P \ where the sum runs over all possible 
configurations of the phases P(x) £ Z3. 

In the context of an effective theory for QCD thermodynamics one con- 
siders the d = 3 model and the spin variables P(x) are interpreted as static 
quark sources at a spatial position x (local Polyakov loops). A vanishing ex- 
pectation value P = V~ 1 ( J £2 X P(x)), which we here will often refer to as the 
magnetization P, signals confinement of quarks, while P / corresponds to 
deconfined quarks. The parameter r is a monotonically increasing function 
of the QCD temperature, k is proportional to the number of flavors and 
decreases with increasing quark mass. The value k = then corresponds to 
pure gluodynamics. [x has the interpretation of the quark chemical potential 
in units of the inverse temperature. 

For vanishing external field, k = (i.e., 77 = fj = 0), the 3-state 
Potts model in three dimensions is known to have a first order transi- 
tion at r = 0.183522(3) [10} I15j . For small external field and vanishing 
chemical potential (i.e., rj = rj = k) the first order transition persists, 
giving rise to a short first order line which ends in a critical endpoint at 
(t,k) = (0.183127(7), 0.00026(3)) [10J. The case of non- vanishing chemical 
potential has been analyzed with techniques based on the Swendsen-Wang 
cluster algorithm using improved estimator techniques in [8] and reweighting 
in [12]. Within the flux representation [7] also local Metropolis updates were 
used [H [TTJ. It has been demonstrated that turning on the chemical poten- 
tial mildens the transition and shifts the critical endpoint towards smaller 
values of k. 

In this article we show that in the flux representation [7] the model is 
accessible also with suitably constructed worm algorithms. We discuss two 
different types of generalized worm algorithms: Type-I: A simple general- 
ization of the conventional Prokof'ev-Svistunov worm [lj with steps where 
monomers are inserted, combined with a random hop to a lattice site with 
another monomer insertion. Type-II: Here the worms are open strings with 
monomers at their endpoints (contrary to the Type-I worms which are closed 
but with the possibility of intermediate hops). 

The two algorithms are compared and evaluated using bulk observables: 
the magnetization P, the internal energy U, the susceptibility xp an d the 
heat capacity C. 



2 



2 Flux representation 



In this section we briefly review the flux representation [7] that we use 
for the two worm algorithms to set the notation and to discuss the flux 
representation of the observables. For the nearest neighbor terms of ([1]) we 
use the ansatz 

e T[P( X )P( X +uy + c.c] = C ^ (P{x)P(x + vT) bX ' V . (3) 

b x , v= -1 

In this expression the term living on the link (x, v) is written as a sum over 
a dimer variable b xv G {— 1,0, +1}. A straightforward calculation gives the 
constants C and B, 



e 2r + 2e- T „ e 2T - e" 

e 2T + 2e" 



C = ^f~ , B = ^- t . (4) 



For the magnetic field terms, which live on a single site x, we use a similar 
ansatz, 

+i 

e nP(x) + v P(xy = J- M Sm P(x) s * , (5) 

Sx = — 1 

which expresses the lhs. as a sum over a monomer variable s x G { — 1, 0, +1} 
attached to the site x. Again one easily works out the monomer weights M s 
for s = —1, 0, +1 and obtains 



M s = - 
3 



(6) 



The weights M s turn out to be non-negative for all values of r\ and r]. 

With the expressions (|3|) and ([5]) the partition sum assumes the form (V 
again denotes the total number of lattice points) 



{6} { S } V x 



,1/ ^-y^JT^s 



(7) 



The first two sums in (JT]) run over all configurations of the dimer variables 
b X)fl and of the monomer variables s x . We have introduced the weight factor 



w[b, S ] = (n^ 1 ) (jiM Sx j . 



(8) 



The summation over the phases P(x) can now be performed in closed form. 
Using that the sum over the roots of unity vanishes, one finds (n is an 
integer), 

Y,P n = 3T(n) with T(n) = 5 nmod ^ , (9) 
p 

where we introduced the triality function T{n) which equals 1 if n is a 
multiple of 3 and vanishes otherwise. The partition sum in its final form 
reads 

z = (3cYEE w MlI r (E[ & ^- k-^l+O • ( 10 ) 

{6} {s} x \ v J 

In ()10p the Potts model is represented in terms of the dimer variables b x v £ 
{ — 1,0, +1} which live on the links (x,u) of the lattice, and the monomer 
variables s x 6 {— 1,0, +1} attached to the sites x. Each configuration of 
dimers and monomers comes with a real and non- negative weight W[b, s] 
as given in (|8|) and thus the complex phase problem is solved. The weight 
consists of a factor B for every non-vanishing dimer, and a monomer factor 
M Sx for each site x, according to the monomer variable s x at that site. 

The configurations of dimers and monomers must obey the constraint 
given by the second factor in ()10p . This factor forces the total flux of dimers 
and monomers at each site x to be a multiple of 3. 

In Fig. [T] we show (for the 2-dimensional case) a few admissible flux 
configurations at individual sites: dimers are represented by arrows in the 
horizontal plane and monomers by arrows on the vertical axis. If an arrow 
points towards a site (away from the site) it contributes with +1 (—1) to 
the total flux at that site. Trivial dimers, b xv = 0, and trivial monomers, 
s x = 0, do not contribute to the constraint and thus need not be represented 
in the illustration of admissible fluxes. 




Figure 1: Examples of admissible fluxes at a site x for the case of d = 2 
dimensions. See the text for the description of the representation. 

In our tests of the worm algorithms we use the bulk observables: internal 
energy U, the magnetization P, i.e., the expectation value of the averaged 
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spin variable P = V~ l J2 x (P( x ))i the heat capacity C and the susceptibility 
Xp- These quantities may be obtained as suitable derivatives of InZ. For 
the evaluation of these derivatives it is useful to note the following relations 
for derivatives of the monomer weights M s : 

dM +1 _ <9M_! _ dM _ dM +l _ dM^ _ dM _ 

drj dfj ' drj dfj 1 ' drj dfj +1 

(11) 

Exploiting these relations one obtains for the expectation value of the spa- 
tially averaged spin the following representation in terms of monomers 

1 d ,.„ ry 1 / <, M o , c M -i c M +i 



P = In Z = — ( 5+i — - + S — — + 5_i — ) , (12) 

V dr] V \ M+i u M M_i / ' v ' 

where <S+i denotes the number of sites x where the monomer variables have 
the value s x = +1, etcetera. 

In a similar way one obtains for the internal energy 



U EE (H) = 



d d _ d 
dr ^ dr/ ^ drj. 



InZ (13) 



= -2dVrB - (BB + S +1 M +1 + S M + 5_i M_i 

with ,6 denoting the number of links (x, z/) where the dimer variables are 
non-zero, b X)V / 0. The new numerical constants in this expression are 

B = 5 9T „ , M +1 - 



e 3r + 1 _ e -3r ' +^ M+1 

~ qM-i+rjM+i ~ r/M +1 +r?M 

M ° = ' M - X = • (M) 

Using equivalent steps also the heat capacity C and the susceptibility \P 
can be expressed in terms of the monomer and dimer numbers and their 
fluctuations. For the evaluation of spin correlators one may use locally 
varying weights r\ and rj as generating functionals. The spin correlators are 
then related to correlators of monomers. 
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3 Worm algorithms 



In a worm algorithm [T], the constraint of flux conservation is temporarily 
relaxed during a Monte Carlo move. One unit of flux is inserted at some 
lattice site and then moved through the lattice by local steps, until the 
insertion is healed again, either by coming back to the original site, or by 
a change of the local monomer number. The individual local steps are 
performed with Metropolis update probabilities. Each complete worm makes 
up one update step in the model with constraint ([9]) which is to be simulated. 

As outlined in the introduction, we use two different types of worm 
algorithms. Type-I: The worms are closed, and in addition to the usual 
local worm propagation, there are random hops with insertion of monomers. 
Type-II: The worms are open strings on the lattice with monomers at the 
endpoints. 

3.1 Type-I: Closed worm algorithm 

Each worm is generated using four different moves. We illustrate them in 
Fig. [2] for the simplest case of an initially empty lattice, i.e., all dimer- 
and monomer variables are set to in the beginning. The worm starts at 
a random position of the lattice (1). It may decide to insert dimers and 
move to the neighboring site (2) but also monomers can be inserted (3). 
The insertion of a monomer is followed by a sequence of attempts to insert 
another monomer at a randomly chosen site, which will succeed (3) with a 
probability depending on the weights M s . These steps are continued until 
the worm closes (4). The acceptance of each step is governed by a Metropolis 
decision. 

The subsequent pseudo-code describes the algorithm for the general case 
of d dimensions (for the numerical tests discussed later we set d = 3) . By xq 




Figure 2: Example of a Type-I closed worm on an initially empty 2-D lattice. 
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we denote the starting point of the worm and x is the current position of its 
head. worm_sign is a variable that determines whether the worm raises or 
lowers the flux, insert-monomer -flag is a logical variable controlling the 
status of the worm. By x © y we denote addition modulo 3 which is the 
usual addition operation except in the cases x = y = 1 and x = y = — 1 
where we define +1 © +1 = —1 and —1 © —1 = +1. Finally randQ is 
a random number generator for uniformly distributed real numbers in the 
interval [0, 1). 

Pseudocode for Type-I closed worms: 
select a lattice site xo randomly 

X i Xq 

select wormsign £ { — 1, +1} randomly 
insert -monomer -flag < — false 

repeat until worm is complete: 

if insert-monomer-flag is true then 
select a lattice site x' randomly 

s < — s x i © {—wormsign) 
if randO < M~JM SJ then 
s x i < — s 

X i X 1 

insert-monomer-flag < — false 
end if 
else 

select i/G{0, ±1,±2, ... ± d} randomly 
if v = then 

s < — s x © wormsign 

if randO < M~/M Sx then 

S x i S 

insert-monomer -flag < — true 
end if 
else 

b < — b x ^ v © [sgn(y) x wormsign) 
if randO <: 5l 6 M^."l then 
b x , v < — b 

X i X + V 

end if 
end if 
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end if 

if x = xq and insert jmonomer_ flag is false then 

worm is complete 
endif 

end repeat until worm is complete 

It is straightforward to show detailed balance using the Boltzmann weight 
W[b, s] given in (JSj) , and that the algorithm is ergodic. 

3.2 Type-II: Open worm algorithm 

We also employed a second variant of the algorithm which we refer to as 
Type-II open worms. In this case, each worm is generated using three moves 
which we again illustrate for an initially empty lattice (Fig. [3]): The worm 
tries to start at a random position of the lattice by inserting a monomer 
(1). It then may decide to insert dimers or monomers. If a dimer is inserted 
the worm moves to the neighboring site (2). If a monomer is inserted, the 
worm is complete and terminates (3). Again the acceptance of the moves is 
decided with a Metropolis step. 

We stress that the Type-II open worm will not work in the case of vanish- 
ing external magnetic field, i.e., for k = 0. In that case no monomers appear 
and the Type-II worm cannot start. For small k the Type-II open worm will 
have a very small probability for starting, since the monomer weights are 
linear in k (for small values of k) . For larger values of k there is no starting 
problem of the Type-II worms (see also Section 4.2). For k ^ it is straight- 
forward to show detailed balance and ergodicity for Type-II worms, and we 
found that they are even simpler to implement than the closed worms. 




• • • • • • 

Figure 3: Example of a Type-II open worm on an initially empty 2- 
dimensional lattice. See the text for its discussion. 
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Figure 4: Phase diagram [14J as obtained from the maxima of the suscepti- 
bility of P. We show results for three values of k. The dashed curves at the 
bottom are the results of a perturbative expansion in r and the horizontal 
line marks the critical value of r for the k = case. 

4 Numerical assessment of the worms 

In this section we show the results of our tests and performance analysis 
for the two generalized worm algorithms. For a detailed discussion of the 
physics results obtained in the simulations of the generalized Potts model 
we refer the reader to |14| . Here we only include a figure of the phase 
diagram to allow the reader to locate the parameter values of the simulations 
discussed in this section in the phase diagram. In Fig. [4] we show the phase 
diagram in the r-fj, plane obtained from the position of the maxima of xp 
for three values of k. The dashed curves at the bottom are the results of a 
perturbative expansion in r and the horizontal line marks the critical value 
of t for the k = case. To the left of the phase lines the Polyakov loop 
is small (confinement), while it is close to 1 to the right of the phase lines 
(deconfined phase). The transitions are of a crossover nature. 

4.1 Validity of the algorithms and the flux-based observables 

As a first test we compare the results from the two worm algorithms with the 
results from an exact evaluation on a small 2x2x3 lattice using the original 
spin formulation of the model. We consider the observables U, C, P and xp 
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Figure 5: Comparison of the exact result (full curves) on a 2 x 2 x 3 lattice 
with the two worm algorithms (symbols). We show the results for U,C,P 
and xp as a function of r using k = 0.001 and ji = 6.0. Note that the 
data points from the two algorithms agree almost perfectly such that the 
corresponding symbols essentially fall on top of each other. 

normalized by the volume. In Fig. [5] we show these observables as a function 
of r for k = 0.001 and /x = 6.0. The curves represent the results from the 
exact evaluation, while for the data from the worm algorithms symbols are 
used. It is obvious, that the numbers from the exact evaluation and from 
the worm algorithms agree very well. This not only establishes the validity 
of the worm algorithms, but also of the mapping of the physical observables 
onto the flux representation. The test on the small volume was repeated for 
different values of the parameters r, k and fi and we always found excellent 
agreement of the worm results with the exact evaluation of the observables. 

In addition to the comparison to exact results on small volumes we per- 
formed also a comparison with a conventional Metropolis calculation in the 
spin representation on larger volumes. This comparison is of course limited 
to vanishing chemical potential, /i = 0, where the complex phase problem is 
absent also in the spin representation and the conventional Metropolis ap- 
proach is possible. In Fig. [6] we compare the results for the internal energy 
U and the heat capacity C as a function of r calculated with the three al- 
gorithms: Conventional Metropolis (red circles), closed worm (blue upward 
pointing triangles) and open worm (green downward pointing triangles) for 
k = 0.01 and \i = on a 16 3 lattice for values of r near the crossover. The 
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Figure 6: Internal energy U (left) and heat capacity C (right) for k = 
0.01 and /i = on a 16 3 lattice at values of r near the crossover. We 
compare the results from the conventional Metropolis simulation in the spin 
representation (red circles), the Type-I closed worm algorithm (blue upward 
pointing triangles) and the Type-II open worms (green downward pointing 
triangles). 



Metropolis simulation uses 10 6 sweeps of local Metropolis steps for equilibra- 
tion and 10 6 measurements separated by 20 sweeps for decorrelation. For 
both worm algorithms we used 10 6 worms for equilibration and 10 6 mea- 
surements. For the Type-I closed worm we separated two measurements by 
Ndeco = 20 worms for decorrelation, while for the Type-II open worms 5000 
worms were used for decorrelation. The results from the three algorithms 
fall on top of each other. The same tests were conducted at other values 
of the parameters and also for P and xp- Furthermore, for [i > where 
simulations with the conventional spin-based Metropolis algorithm are not 
possible, we performed a comparison among the two worm algorithms and 
also in this case found excellent agreement. We conclude that the two worm 
algorithms are valid, and that the representation of the observables in terms 
of fluxes was implemented correctly. 

4.2 Time series 

In this subsection we now study performance aspects of the two worm algo- 
rithms. We begin our assessment with analyzing the time series of measure- 
ments of the total number B of non-zero dimer variables and the number 
S of non-zero monomers. As shown at the end of Section 2, the physical 
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Figure 7: Time series for B/V and S/V at r = 0.1782, k = 0.005, and /i = 
on a 8 3 lattice, generated with Type-I worms. The Monte Carlo time on the 
horizontal axis is measured in units of A^ eC o = 20 worms. 



observables can be constructed from B and S and their moments. In Fig. [7] 
we show time series from Type-I worms for B/V in the top plot and for S/V 
(bottom) for the parameters r = 0.1782, k = 0.005, \x = 0.0 and a volume 
of 8 3 (time series for Type-II worms behave similarly). Between two mea- 
surements we used Nd eco = 20 worms for decor relation, so the Metropolis 
time on the horizontal axis is measured in units of 20 worms. We do not 
observe structures on scales larger than several hundred worms and as a 
first impression conclude that the Type-I worms lead to reasonably quick 
sampling. The sampling becomes better for larger values of k or /i. 

There are, however, regions in the space of couplings where one finds 
that the flux representation and worms are not the optimal approach. In 
general this is the case for very small values of k (typically for k of O(10 -4 )), 
where the monomer weights are small (in leading order they are linear in n) , 
and thus it is hard to generate monomers. This is illustrated in Fig. [HI where 
we show the same time series as in Fig. [7j but now at k = 0.0002, /i = 
and t = 0.183, again with Nd eC o = 20 worms for decorrelation. Although 
the time series for B/V still looks reasonable, the time series for S/V shows 
that monomers are generated only very rarely and very high statistics would 
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Figure 8: Time series for B/V and S/V at k = 0.0002, /u = and r = 0.183 
on a 8 3 lattice, generated with Type-I worms. The Monte Carlo time on the 
horizontal axis is measured in units of N^eco = 20 worms. 



be needed to obtain reliable results. We stress that this is a property of the 
flux representation: As mentioned, the monomer weights are very small 
at small k, so monomers are indeed expected to be very rare and the flux 
representation is not the optimal choice for very small k. This difficulty could 
be overcome by reweighting the sector of phase space with small monomer 
numbers to obtain a higher weight and thus higher frequency in the Monte 
Carlo. They could then be weighted back in order to compute expectation 
values. In the present study we did not explore such a reweighting strategy. 
We also note that for small k the complex phase problem is very mild and 
reweighting techniques were used to successfully simulate the system with 
the conventional spin representation [U |T2] in that region of the parameter 
space also at small fi. 



4.3 Characteristic quantities 

Let us now come to a quantitative comparison of the Type-I closed worms 
and Type-II open worms. For a sensible comparison one has to consider 
the actual changes done by the worms: A Type-II open worm inserts two 



13 



monomers at its endpoints and connects it with dimers (one segment), while 
a Type-I closed worm may consist of several segments of dimer strings con- 
nected by large hops between monomers (it may also consist of only dimer 
hops without large hops, i.e., only one closed segment). We first look at some 
characteristic quantities in order to describe the behavior of the worms in 
different regions of parameter space. By n os we denote the average number 
of open segments in Type-I worms, i.e., segments with monomer hops at 
their endpoints (if a Type-I worm closes without monomer insertions tl os 
is zero). By D we denote the average number of dimer steps in a Type-I 
or Type-II worm. Another interesting quantity for the comparison of the 
two worms is the ratio r of the number of successful worm starts divided by 
the number of all start attempts. Depending on the parameters, here one 
expects a drastic difference between the worms, since the open worms can 
start only with the insertion of a monomer. Finally, for the comparison of 
the computational cost of the two algorithms we define the cost ratio c s for 
the Type-I and Type-II worms as the ratio of the total number of attempted 
steps (dimer, monomer and start attempts) to the number of accepted steps. 

In Table [U we collect the data for Type-I closed worms, and in Table 
[2] for Type-II worms. The analysis was done for five different sets of pa- 
rameters (first column) and for four volumes (second column). Let us begin 
with the comparison of the ratio r (number of successful worm starts di- 
vided by the number of all start attempts). For Type-I closed worms the 
starting probability is similar for all parameter sets (ranging between 0.129 
and 0.209). Since here the chance for starting with a dimer insertion is 6/7, 
the ratio goes down slightly when dimers are more costly, i.e., when t is 
small, or when the system is monomer dominated (k and/or fi large). For 
the Type-II open worms, r shows a much larger variation: Here the worm 
must start with the insertion of a monomer, and for some parameter sets 
(Sets B, D and E) these are scarce, such that for the Type-II worms many 
start attempts fail. On the other hand, for parameters where monomers are 
abundant (Sets A and C), we observe a very high rate of successful starts. 

For Type-I closed worms the abundance of monomers is also reflected in 
the average number n os of open segments: It is small for the same parameter 
sets (Sets B, D and E) where the starting ratio r of the Type-II worms is 
small, i.e., for regions of parameter space with few monomers. Concerning 
the dependence on the volume V, we find that n os grows linearly with the 
volume, i.e., n os = cV perfectly describes the data, where c is a constant 
that depends on the parameters r, k and [i. 

The behavior of the average number D of dimer steps is interesting. 
For Type-II open worms we find that D varies only very little with the 
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Table 1: Characteristic quantities for Type-I closed worms (see the text for 
their definitions). All the simulations use 10 6 worms for equilibration and 
10 7 measurements separated by N^ eco = 5 worms for decorrelation. The 
error is usually smaller than the last digit we show. 

volume and is correlated with the abundance of monomers: Since monomers 
are needed to terminate Type-II worms, parameter values where monomers 
are scarce (Sets B, D and E) will lead to longer worms, i.e., increase D. 
For Type-I closed worms we find that D essentially grows linearly with the 
volume for Sets A, C and E. For the sets B and D (those with very low density 
of monomers) we find that D grows slower than linear. This is correlated 
with a very small number n os of open segments for these ensembles. 

The cost ratio c s essentially reflects the behavior of r, n os and D: When- 
ever monomers are sufficiently abundant, the values of c s are independent 
of the volume and are roughly equal for Type-I and Type-II worms. Dis- 
crepancies are observed only for the sets where monomers are scarce (Sets 



15 



Parameters 


T/ 
K 


T 


u 


c s 


set a: 





n orq 

U.ZDO 


z.yo 


Q 7Q 

o. to 


T — U.1UU, 





U.ZDO 


z.yo 


Q 7Q 
O. (O 


„ _ n nni 

rt — U.UUlj 


1 fi 3 

-LU 


D 9fi8 

U.iUO 




o. / o 


., C Q 


qo3 

OZ 


n 9fi8 
U.ZDO 


9 d 1 } 

z.yo 


'i 78 
0. lO 


set d. 


D 


n nm a 

U.UU14 


1 ns 

1U0 


Q AA 
y .44 


i — nisi 

T — U.101, 


8 3 
o 


n nm k 

U.UU10 


1 

IDo 


7 17 
I.LI 


K — n nni 


1 fi 3 

-LU 


D 001 fi 

u. uu ±u 


1 Q^i 

-L tju 


fi 49 


/i — U.O 


oo3 
OZ 


n nm ^ 

U.UU10 


1 70. 

i /y 


R 77 
D. / 1 


set 


D 


U.DZO 


n /i /i k 

U.44D 


o.OD 


n no^ 

T — U.Uzo, 


«3 
o 


n 

U.OZO 


C\ A AK 
U.44D 


O.OO 


k — 005 


16 3 


0.625 


0.445 


3.85 


(i = 5.2 


32 3 


0.626 


0.445 


3.85 


Set D: 


6 3 


0.0059 


30.1 


8.57 


r = 0.170, 


8 3 


0.0058 


29.2 


8.77 


k = 0.005, 


16 3 


0.0058 


28.3 


8.98 


H = 0.2 


32 3 


0.0058 


28.3 


8.98 


Set E: 


6 3 


0.0482 


12.5 


4.70 


r = 0.150, 


8 3 


0.0482 


12.5 


4.70 


k = 0.010, 


16 3 


0.0482 


12.5 


4.70 


/i = 2.0 


32 3 


0.0482 


12.5 


4.70 



Table 2: Characteristic quantities for Type-II open worms (see the text for 
their definition). All the simulations use 10 6 worms for equilibration and 
10 7 measurements separated by N^ eco = 500 worms for decorrelation. The 
error is usually smaller than the last digit we show. 

B and D): In these cases the closed Type-I worms have the advantage that 
they can run without any monomer insertions at all, thus improving their 
cost ratio. This advantage is, however, lost when the volume is larger and 
longer dimer chains emerge. Then from time to time a monomer insertion 
takes place and also Type-I worms then spend quite some time with trying 
to insert another monomer necessary for continuation. Thus for larger vol- 
umes the cost ratios c s of Type-I and Type-II worms roughly agree also for 
Sets B and D. From the assessment of r, n os , D and c s we conclude that for 
these characteristic quantities the two algorithms behave similarly, with the 
exception that the closed Type-I worms are more flexible since they work 
also for vanishing external magnetic fields. 
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4.4 Autocorrelation times 

The relevant figure of merit for an algorithm is computational effort for a 
fixed precision of results. We therefore analyze the integrated autocorrela- 
tion time r^ t of the bulk observables X = U,C, P and \p- Since the sizes 
of worms vary drastically between the different cases, we need to normal- 
ize the autocorrelation times. For this purpose we define one "sweep" as 
ro = 3V/D worms, i.e., the average number of worms needed to visit every 
link of the lattice as the (customary) unit for the integrated autocorrelation 
times r^ t . Obviously one may express ro also in terms of measurements as 
ro = 3V/(D Ndeco) measurements. In Tables 3 and 4 we give the autocorre- 
lation times in units of To- In order to obtain a measure of computational 
effort, the results are multiplied by the cost ratio c s , in other words we 
show T = c s T meas /TQ, where T meas simply is the unnormalized autocorrela- 
tion time in units of measurements. The statistical errors of autocorrelation 
times were estimated by a jackknife procedure and are between about 5 and 
10 percent depending on the parameter values. This is sufficient for the 
subsequent comparison of the two algorithms. 

The autocorrelation times are almost independent of volume for most 
parameter sets, reflecting the crossover or noncritical nature of correlations. 
Even for set B, which is very close to the critical value of r at k = 0, 
autocorrelation times of Type-I are small and increase only very moderately 
with volume. Self averaging within a large volume can lead to the decreasing 
times (at high scale) for the internal energy U for set B and Type-II. 

Comparing the autocorrelation times (Tables 3 and 4) shows that Type- 
II open worms outperform the Type-I closed worms (i.e., shorter autocorre- 
lation times) for those sets where monomers are abundant (Sets A and C), 
while the Type-I worms perform better for Sets B, D, E, where monomers 
are scarce. The reason is that Type-II worms always start with the insertion 
of a monomer, while Type-I worms insert monomers only with a probability 
1/7. Thus Type-II worms may have a larger portion of monomer insertions 
compared to dimer steps, and thus do many monomer changes which is ad- 
vantageous for Sets A and C. Vice versa, for the sets with few monomers, 
the Type-I worms can do many dimers steps without needing monomer in- 
sertions and thus perform better in this situation. We remark that some of 
the disadvantage of the Type-I worms for parameter sets where monomers 
are abundant could be ameliated by changing the probability for offering a 
monomer change from 1/7 to some larger value. Remarkably, the relative 
performance of the two worm types is roughly proportional to r/c s , i.e. to 
the acceptance rates for attempted steps. 
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Table 3: Type-I closed worms: Autocorrelation times of the bulk observables 
at different parameters. We show the autocorrelation times f in units of tq 
and multiplied by the cost ratio c s as discussed in the text. 

5 Summary and discussion 

In this article we presented in detail two worm algorithms that were used in 
a recent QCD-inspired study of the 3-state Potts model with external field 
and chemical potential in three dimensions [13] and analyzed their perfor- 
mance. The algorithms are based on a flux representation of the model which 
can be obtained using high temperature expansion techniques. The nearest 
neighbor term leads to the conventional dimer-based closed contours, while 
the magnetic terms allow violation of the constraint through the insertion 
of monomers (similar structures will appear whenever one considers spin 
systems that are coupled to an external field). Furthermore, the chemical 
potential gives a different weight to monomers and antimonomers. Thus the 
model provides a rather general testbed for various algorithmic ideas based 
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Table 4: Type-II open worms: Autocorrelation times of the bulk observables 
at different parameters. We show the autocorrelation times f in units of tq 
and multiplied by the cost ratio c s as discussed in the text. For cases marked 
with ">", only a lower bound could be determined. 

on the worm concept PQ. 

We explored two types of generalized worm algorithms that differ in 
their treatment of monomers. The Type-I closed worms start with dimers 
or monomers and allow for dimer insertions or the insertion of monomers 
which are then followed by a random hop to another position where a second 
monomer is inserted. These closed worms finish when the starting position 
is reached again. On the other hand the Type-II open worms create chains 
of dimers that start and end with the insertion of monomers. It is important 
to note that the open worms are restricted to simulations with non-vanishing 
external magnetic field. 

It was carefully checked that the two worm algorithms produce correct 
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results, by comparison with an exact evaluation on small volumes, a cross- 
check with a standard Metropolis simulation (at vanishing chemical poten- 
tial) , and a systematic comparison of the results from the two algorithms at 
various parameter values. 

In order to evaluate the performance of the two algorithms we first com- 
pared various characteristic quantities, namely the starting probability r, 
the number of open segments n os in Type-I worms, the average number of 
dimer steps D of the worms, and their cost ratio c s . We then presented 
the results of an analysis of the autocorrelation times, suitably scaled to 
reflect computational effort. The overall assessment shows that the Type-I 
closed worms are more flexible: They can be used also at vanishing external 
field and in general were found to perform better for parameter sets where 
monomers are scarce. Type-II open worms perform considerably better in 
monomer dominated regions of the parameter space. However, we expect 
that changing the ratio of monomer to dimer steps in Type-I closed worms 
would improve the latter for situations where monomers are abundant. 

The authors expect that generalizations of the worm concept will see 
a lot of attention in the future: Matter fields (bosonic and fermionic) give 
rise to closed loops of flux on a lattice, and worm algorithms are a natural 
approach to update such a system. Of course in more realistic models the 
structure of the fluxes is more involved (see, e.g., the generalized effective 
Polyakov loop model [16]) and new strategies need to be found. We expect 
that the current paper is useful for such further developments. 
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